Recall the test problem

$$v'''(t) + v'(t)v(t) - \frac{\beta_1 + \beta_2 + \beta_3}{3}v'(t) = 0,$$

where $\beta_1 < \beta_2 < \beta_3$. It follows that

$$v(t) = \beta_2 + \left(\beta_3 - \beta_2\right)\text{cn}^2\left(\sqrt{\frac{\beta_3 - \beta_1}{12}}t, \sqrt{\frac{\beta_3 - \beta_2}{\beta_3 - \beta_1}}\right)$$

is a solution where $\text{cn}(x, k)$ is the Jacobi cosine function and $k$ is the elliptic modulus. Some notations use $\text{cn}(x, m)$ where $m = k^2$. The corresponding initial conditions are

$$v(0) = \beta_3, \quad v'(0) = 0, \quad v''(0) = -\frac{\left(\beta_3 - \beta_1\right)\left(\beta_3 - \beta_2\right)}{6}.$$

Let's write the equation as a system and compute the Jacobian. For $\beta_1 = 0$, $\beta_2 = 1$, and $\beta_3 = 10$, based on an analysis of the Jacobian, we will see if we can suggest methods to solve the problem.

Letting $u_1 = v(t)$, $u_2 = v'(t)$, and $u_3 = v''(t)$, we can express the nonlinear system as $u' = f(u)$ where

$$f(u) = \begin{bmatrix} u_2 \\ u_3 \\ \frac{\beta_1 + \beta_2 + \beta_3}{3}u_2 - u_1u_2 \end{bmatrix}.$$

The Jacobian, $f'(u)$, is then given by

$$\begin{bmatrix} 0 & 1 & 0 \\ 0 & 0 & 1 \\ -u_2 & \frac{\beta_1 + \beta_2 + \beta_3}{3} - u_1 & 0 \end{bmatrix}.$$

The range of possible eigenvalues for different values of $t$ are plotted below.

Due to the location of the eigenvalues it is not possible (see HW3!) to make it so that all the eigenvalues lie within the stability region for any method. So, it is about maximizing. With something like trapezoid, the eigenvalues are in the stability region half of the time. With backward Euler, this is true of y2, but the others are within the stability region more than half of the time.

Something like leapfrog would be a bad choice here.

The method one should choose for this problem is maybe less clear. But similar considerations as those given for the previous choice of $\beta_i$'s still are good justification. But for both problems, Runge--Kutta 4 is actually probably the best choice -- one-step, accurate, and for sufficiently small step size the eigenvalues will be within the stability region approximately half the time.